1. Field of the Invention
The present invention relates to an image enhancement technique. More particularly, the present invention relates to an image enhancement method and system that are capable of improving image quality by performing gamma corrections on global and local illuminations and reflectance estimated from an input image, respectively.
2. Description of the Related Art
With the widespread popularity of digital cameras (including digital camera-equipped mobile phones), taking images has become as commonplace as making a telephone call. However, most images, such as those obtained with a portable digital camera, are obtained in an irregular illumination environment such as fluorescent lamp, sunlight, and street lamp, rather than in an artificially prepared regular illumination environment, e.g. in a studio. Such an image is likely to be deteriorated with an excessively shaded or bright region. Accordingly, various techniques have been proposed for improving the quality of the deteriorated image.
Well known techniques for improving image quality include intensity transformation, histogram modeling, and homomorphic filtering, retinex and Kimmel methods based on an image formation model.
The intensity transformation technique is disclosed by R. C. Gonzalez and R. E. Woods in “Digital Image Processing,” Reading, Mass., Addison-Wesley, 1992 and W. K. Pratt, “Digital Image Processing,” 2nd ed. New York, Wiley, 1991.
The histogram modeling technique is disclosed by A. K. Jain in “Fundamentals of Digital Image Processing,” Enblewood Cliffs, N.J., Prentice-Hall, 1989.
The homomorphic filtering method based on the image formation model is disclosed by J. S. Lim in “Two-Dimensional Signal and Image Processing” Englewood Cliffs, N.J., Prentice-Hall, 1990.
The retinex method is disclosed by D. J. Jobson, Z. Rahman, and G. A. Woodell in “Properties and Performance of a Center/Surround Retinex,” IEEE Trans. Image Process., vol. 6, no. 3, pp. 451-462, March 1997, by M. Ogata, T. Tsuchiya, T. Kubozono, and K. Ueda in “Dynamic Range Compression Based on Illumination Compensation,” IEEE Trans. Consumer Electron., vol. 47, no. 3, pp. 548-558, August 2001, by R. Kimmel, M. Elad, D. Shaked, R. Keshet, and I. Sobel in “A Variational Framework for Retinex,” Int. J. Comput. Vis., vol. 52, no. 1, pp. 7-23, January 2003, and by D. J. Jobson, Z. Rahman, and G. A. Woodell in “A Multi-Scale Retinex for Bridging the Gap Between Color Images and the Human Observation of Scenes,” IEEE Trans. Image Process., vol. 6, no. 7, pp. 965-976, July 1997.
The Kimmel method is disclosed in U.S. Pat. No. 6,941,028 and by R. Kimmel, M. Elad. D. Shaked, R. Keshet, and I. Sobel in “A Variational Framework for Retinex,” Int. J. Comput. Vis., vol. 52, no. 1, pp. 7-23, January 2003.
All of these methods improve the qualities of images by decreasing a dynamic range of an input image or increasing a contrast of the image.
The intensity transformation method, histogram modeling method, homomorphic filtering method, retinex method and Kimmel method will each be briefly described below. The intensity transformation method improves the quality of an input image by transforming brightness values of the image using a linear function, log function or power function (R. C. Gonzalez and R. E. Woods, “Digital Image Processing,” Reading, Mass., Addison-Wesley, 1992 and W. K. Pratt, “Digital Image Processing,” 2nd ed. New York, Wiley, 1991). In the case of using a linear function as the transformation function, it is called contrast stretching or gain/offset correction method. In the case of using a power function, it is called gamma correction method. These intensity transformation methods are known to be easily implemented. However, since only one transformation function is used for an entire image, it is difficult, when there exist several regions in the image of which contrasts are different from each other, to increase the contrasts at the same time.
The histogram modeling method obtains an improved image by transforming the input image into an image which has a required form of histogram. In the case that the required form of histogram is uniform, it is called a histogram equalization method. Although the histogram equalization method can be easily implemented as a brightness transformation, it is difficult to increase the contrast of different areas in an image at the same time. For example, the contrast of a simple background distributed in a small area decreases, while the contrast of a simple background distributed in a large area increases. Accordingly, the processed image is unnaturally expressed rather than improved.
In order to overcome this shortcoming, a modified histogram equalization method has been proposed in which an input image is divided into a plurality of blocks and the histogram equalization is applied for each block. Such an algorithm is disclosed by R. C. Gonzalez and R. E. Woods in “Digital Image Processing,” Reading, Mass., Addison-Wesley, 1992; by S. M. Pizer, E. P. Amburn, J. D. Austin, R. Cromartie, A. Geselowitz, T. Greer, B. T. H. Romeny, J. B. Zimmerman, and K. Zuiderveld in “Adaptive Histogram Equalization and its Variations,” Comput. Vis., Graph., Image Process., vol. 39, no. 3, pp. 355-368, September 1987 and by J. Y. Kim, L. S. Kim; and S. H. Hwang in “An Advanced Contrast Enhancement Using Partially Overlapped Sub-block Histogram Equalization,” IEEE Trans. Circuits Syst. Video Technol., vol. 11, no. 4, pp. 475-484, April 2001.
The R. C. Gonzalez and R. E. Woods, “Digital Image Processing” Reading, Mass., Addison-Wesley, 1992, discloses a method in which the input image is divided into blocks of pixels; the S. M. Pizer, E. P. Amburn, J. D. Austin, R. Cromartie, A. Geselowitz, T. Greer, B. T. H. Romeny, J. B. Zimmerman, and K. Zuiderveld, “Adaptive Histogram Equalization and its Variations,” Comput. Vis., Graph., Image Process., vol. 39, no. 3, pp. 355-368, September 1987, discloses a method in which the blocks are divided without overlapping with each other; the J. Y. Kim, L. S. Kim, and S. H. Hwang, “An Advanced Contrast Enhancement Using Partially Overlapped Sub-block Histogram Equalization,” IEEE Trans. Circuits Syst. Video Technol., vol. 11, no. 4, pp. 475-484, April 2001 discloses a histogram equalization method with partially overlapped blocks; and the V. Buzuloiu, M. Ciuc, R. M. Rangayyan, and C. Vertan, “Adaptive-neighborhood Histogram Equalization of Color Images,” J. Electron. Imag., vol. 10, no. 2, pp. 445-459, April 2001 discloses a histogram equalization method in which the image is divided into blocks in consideration of the variations of the size and form of the input image. Also, the V. Buzuloiu, M. Ciuc, R. M. Rangayyan, and C. Vertan, “Adaptive-neighborhood Histogram Equalization of Color Images,” J. Electron. Imag., vol. 10, no. 2, pp. 445-459, April 2001 discloses a generalized histogram equalization method.
The homomorphic filtering method applies a logarithmic function to the input image and then applies a linear low-pass filter (LPF) and a linear high-pass filter (HPF) to the log signal of the input image. From the respective output, the log signal of the illumination and the log signal of the reflectance are estimated. The estimated log signal of the illumination is multiplied by a value less than 1 for reducing the dynamic range and the estimated log signal of reflectance is multiplied by a value greater than 1 for increasing the contrast. Finally, the two signals are summed at an adder and an exponential function is applied to the sum so as to return to the spatial domain from the log domain. In this method, by multiplying a value greater than 1 to the log signal of the reflectance, the contrast of a bright region in the output image may increase in relation to the contrast of a dark region. So, the characteristic of the human visual system (HVS) is not well reflected. The nonlinear log operation generates harmonics in the frequencies of the input image. This results in variation of the frequency spectrum so as to cause unreliability of estimation on the log signals of the illumination and the reflectance.
Meanwhile, the retinex method based on the retinex theory proposed by E. Land in “An Alternative Technique for the Computation of the Designator in the Retinex Theory of Color Vision,” Proc. Nat. Acad. Sci., vol. 83, pp. 3078-3080, May 1986, assumes that the HVS recognizes colors without affection of the illumination.
The retinex method applies a linear LPF for estimating the illumination from the input image and then the reflectance is estimated by removing the estimated illumination from the input image. Only the reflectance is next emphasized. The retinex method proposed by D. J. Jobson, Z. Rahman, and G. A. Woodell in “Properties and Performance of a Center/surround Retinex,” IEEE Trans. Image Process., vol. 6, no. 3, pp. 451-462, March 1997 estimates the illumination using a Gaussian-form linear LPF and subtracts the log signal of the estimated illumination from the log signal of the input image so as to estimate the log signal of the reflectance. Without returning to the spatial domain, the reflectance signal estimated in the log domain is used as the output image. Finally, the brightness range of the output image is adjusted by a gain/offset correction in accordance with that of the output device. In the method proposed by Jobson et al., the reflectance of the output image is estimated with a difference between the input image and an illumination estimated using the linear LPF.
Since the output image is given as the estimated reflectance, the output image is characterized in that a dark region becomes darker and a bright region becomes brighter along the edge at which the brightness tends to vary abruptly such that a halo effect appears. However, in a case that the support region of the linear LPF is narrow, the estimated illumination is smoothed in a small area around the edge. So the halo effect appears in a small area around an edge and the local contrast increases. On the other hand, if the support region of the linear LPF is wide, the estimated illumination is smoothed in large area. At this time, the halo effect appears in a large area and the global contrast increases. In order to increase the local contrast and the global contrast simultaneously, a method using the weighted sum of output images by using three linear LPFs having different support regions has been proposed. This method is disclosed by D. J. Jobson, Z. Rahman, and G. A. Woodell in “A Multi-scale Retinex for Bridging the Gap Between Color Images and the Human Observation of Scenes,” IEEE Trans. Image Process., vol. 6, no. 7, pp. 965-976, July 1997. Such a retinex method is called multi-scale retinex (MSR), whereas one that uses only one linear low-pass filter is called single-scale retinex (SSR).
The method proposed by Kimmel is disclosed by R. Kimmel, M. Elad, D. Shaked, R. Keshet, and I. Sobel in “A Variational Framework for Retinex,” Int. J. Comput. Vis., vol. 52, no. 1, pp. 7-23, January 2003. In this method, the illumination is iteratively estimated using quadratic programming (QP) under a constraint that the illumination of the output image is brighter than or equal to that of the input image over the entire region. The reflectance is estimated by inversing the estimated illumination and multiplying it with the input image. In order to compress the dynamic range of the estimated illumination, the gamma correction is applied and then it is multiplied with the estimated reflectance so as to obtain an output image.
The homomorphic filtering method, retinex (SSR, MSR), and Kimmel method will be described hereinafter with reference to FIGS. 1 to 8 in more detail.
In the conventional image formation model, a grayscale image or a component image of a color image f(x,y) can be expressed as a multiplication of illumination l(x,y) and reflectance r(x,y) in equation (1):f(x,y)=l(x,y)·r(x,y)  (1)
where it is assumed that the illumination l(x,y) is slowly changed and its frequency spectrum is mostly distributed at low frequency band. Also, it is assumed that the reflectance r(x,y) is rapidly changed due to the reflection characteristics on the surface of an object and its frequency spectrum is mostly distributed at high frequency band. In the conventional image enhancement technique based on the image creation model of equation 1, at least one of the illumination and reflectance is estimated. The estimated illumination is then processed so as to reduce its dynamic range and the estimated reflectance is processed to increase its contrast. An output image is obtained by combining these processed components. Next, the output image is used as it is or is corrected by gain/offset correction for matching the brightness value range of an output device.
Homomorphic Filtering Method
As shown in FIG. 1, the log signal log f(x,y) of the input image is obtained by applying the logarithmic function at logarithmic function block 101. Based on the image formation model of equation 1, the log f(x,y) is expressed as the sum of the log signal log l(x,y) of the illumination and the log signal log r(x,y) of the reflectance as in equation (2).log f(x,y)=log l(x,y)+log r(x,y)  (2)
The log signal log l(x,y) of the illumination of equation (2) is estimated by applying the LPF 103 to the log signal log f(x,y) of the input image under the assumption that the frequency spectrum is distributed at low frequency band as in equation (3). The log signal log r(x,y) of the reflectance is estimated by applying the HPF 105 as in equation (4) under the assumption that the frequency spectrum of the log signal is of the reflectance is distributed at high frequency band,log {circumflex over (l)}(x,y)=LPF[log f(x,y)]  (3)log {circumflex over (r)}(x,y)=HPF[log f(x,y)]  (4)
where log {circumflex over (l)}(x,y) is the estimated log signal of illumination, log {circumflex over (r)}(x,y) is the estimated log signal of reflectance. The estimated log signal log {circumflex over (l)}(x,y) of the illumination is multiplied by a constant α less than 1 at the multiplier 107 for decreasing the dynamic range, and the estimated log signal log {circumflex over (r)}(x,y) of the reflectance is multiplied by a constant β greater than 1 at the multiplier 108 for increasing the contrast. The output signals of the multipliers 107 and 108 are summed at an adder 109 and then the output signal of an adder 109 is processed by the exponential function block 111 so as to be reconstructed. At this time, getting back to the spatial domain obtains the output image {circumflex over (f)}(x,y). From the comparison of an input image illustrated in FIG. 2A and an output image illustrated in FIG. 2B, homomorphic filtering method improves the quality of the image. However, the homomorphic filtering method requires the nonlinear log operation such that the log signal f(x,y) of the input image generates harmonics of the frequency spectrum. The harmonics of the frequency spectrum make it difficult to estimate the log signals of the illumination and reflectance. Particularly, the output image {circumflex over (f)}(x,y) obtained using the homomorphic filtering method can be expressed by equation (5).{circumflex over (f)}(x,y)=exp(α·log {circumflex over (l)}(x,y)+β·log {circumflex over (r)}(x,y))  (5)
Rearranging the right terms of the equation 5 with respect to the estimated illumination {circumflex over (l)}(x,y) and the estimated reflectance {circumflex over (r)}(x,y) the output image {circumflex over (f)}(x,y) can be expressed in a form of the multiplication of the two estimated components enhanced by the gamma correction as shown in equation 6:{circumflex over (f)}(x,y)=({circumflex over (l)}(x,y))α·({circumflex over (r)}(x,y))β  (6)
where α<1 and β>1. In equation (6), a gamma correction is applied with a gamma factor β greater than 1 to the reflectance in the range of 0˜1, the contrast of the image increases in the bright region relative to the dark region, as shown in FIG. 2B. This is because the homomorphic filtering method does not consider that the HVS is more sensitive to the contrast of the dark region than that of the bright region.
FIG. 2A illustrates an input image having strong edges at the center part and dark and weak edges at the left and right parts. FIG. 2B illustrates the output image obtained by applying the homomorphic filtering to the input image of FIG. 2A. In this case, the log signal of the illumination is estimated by separately applying a 1-dimensional Gaussian LPF with a window of 1×5 size in horizontal and vertical directions to the log signal of the input image. The log signal of the reflectance is estimated by subtracting the log {circumflex over (l)}(x,y) from the log f(x,y). The constants α and β are 0.5 and 1.5, respectively. From the output image of FIG. 2B, it is shown that the bright region around the tower shows strong contrast relative to the dark region around trees.
Retinex Method
Typically, retinex methods are classified into SSR and MSR.
FIG. 3 is a block diagram illustrating the SSR method proposed by D. J. Jobson, Z. Rahman, and G. A. Woodell in “Properties and Performance of a Center/surround Retinex,” IEEE Trans. Image Process., vol. 6, no. 3, pp. 451-462, March 1997.
In the SSR, the illumination is estimated by applying a linear LPF 301 to the input image f(x,y). The input image is also directed to log function block 305. The output image {circumflex over (f)}(x,y) is obtained at adder 307 from a signal taken by removing a log signal log {circumflex over (l)}(x,y) obtained from log function block 303 of the estimated illumination, from the log of the input image obtained from log function block 305 as shown in equation (9):{circumflex over (f)}(x,y)=log f(x,y)−log {circumflex over (l)}(x,y)  (9)
where {circumflex over (l)}(x,y) is the estimated illumination through a linear low-pass filtering by equation 10.{circumflex over (l)}(x,y)=LPF[f(x,y)]  (10)
The output image {circumflex over (f)}(x,y) of equation 9 can be expressed with respect to the estimated illumination {circumflex over (l)}(x,y) by equation (11):{circumflex over (f)}(x,y)=({circumflex over (l)}(x,y))0·log({circumflex over (r)}(x,y))  (11)
where {circumflex over (r)}(x,y) is the reflectance estimated with
            r      ^        ⁡          (              x        ,        y            )        =                    f        ⁡                  (                      x            ,            y                    )                                      l          ^                ⁡                  (                      x            ,            y                    )                      .  
The SSR represented by equation (11) applies a gamma collection with a factor of α=0 to the estimated illumination {circumflex over (l)}(x,y) such that the affection of illumination is excluded unlike in the homomorphic filtering method. However, the logarithmic function, whose transfer characteristic is similar to the gamma correction with a factor of β=0.3, is applied to the reflectance {circumflex over (r)}(x,y) such that the contrast increases excessively at a dark region.
Finally, in order to obtain final output image {tilde over (f)}(x,y) whose brightness value range is adjusted to that of the output device, a gain/offset correction is applied to the output image {circumflex over (f)}(x,y) by gain/offset correction function block 309 as in equation 12:{tilde over (f)}(x,y)=a[{circumflex over (f)}(x,y)+b]  (12)
where a and b are constants used for the gain/offset correction. The values of a and b affect the final output image {tilde over (f)}(x,y). As disclosed by D. J. Jobson, Z. Rahman and G. A. Woodell in “Properties and Performance of a Center/surround Retinex,” IEEE Trans. Image Process., vol. 6, no. 3, pp. 451-462, March 1997, the same values of a and b are used for all images.
FIG. 4 is a diagram illustrating a 1-dimensional signal presenting an edge component of an input image and estimated illumination and reflectance in the SSR method of FIG. 3. In FIG. 4, a linear LPF smoothes the edge component of the input image 401 to estimate the illumination 402 such that the reflectance 403 estimated by a reflectance estimator 405 abruptly decreases or increases around the edge. This causes a halo effect in which the dark region of the input image becomes darker and the bright region of the input image becomes brighter in the output image. The edge region of the smoothed image in the illumination estimated by the linear low-pass filtering varies in size and smoothing degree according to the support region of the filter. The edge region to be smoothed becomes narrow or wide according to the width of the support region of the filter. At this time, the smoothing degree also increases or decrease and the scale of the halo effect increases or decreases. If the support region of the filter is narrow, the halo effect occurs in a small area so the local illumination can be well estimated, resulting in an increase of the local contrast of the image. Conversely, if the support region of the filter is wide, the halo effect occurs in a large area so the global illumination can be well estimated, resulting in an increase of the global contrast of the image. Since the halo effect and the contrast vary according to the support region of the filter, it is required to adjust the support region of the filter. FIGS. 5A to 5C illustrate SSR output images to the input image of FIG. 2A according to the width of the support region of the linear LPF. The constants a and b for the gain/offset correction in equation 12 have values of 192 and −30, respectively, as in D. J. Jobson, Z. Rahman, and G. A. Woodell, “A Multi-scale Retinex for Bridging the Gap Between Color Images and the Human Observation of Scenes,” IEEE Trans. Image Process. vol. 6, no. 7, pp. 965-976, July 1997. If the support region of the filter is set to 15 pixels, the halo effect in the output image of SSR in FIG. 5A appears at small areas around the roof and the tower. It shows that the local contrast of the image increases but the global contrast does not increase well. FIG. 5C illustrates an output image of SSR when the support region of the filter is set to 500 pixels. In this case, the halo effect appears in wide areas around the roof and the tower. This means that the global contrast increases but the local contrast does not increase. FIG. 5B illustrates an output image of SSR when the support region of the filter is set to 150 pixels. In this case, the halo effect appears to an extent between those of FIGS. 5A and 5C. The MSR is proposed by D. J. Jobson, Z. Rahman, and G. A. Woodell in “A Multi-scale Retinex for Bridging the Gap Between Color Images and the Human Observation of Scenes,” IEEE Trans. Image Process., vol. 6, no. 7, pp. 965-976, July 1997. The MSR is a technique for compensating the variational characteristics of the SSR output images according to the support region of the linear LPF for estimating the illumination. In the MSR, an output image {circumflex over (f)}MSR(x,y) is obtained by a weighted sum of the SSR output images taken by using the linear LPFs having different support regions as expressed in equation 13:
                                                        f              ^                        MSR                    ⁡                      (                          x              ,              y                        )                          =                              ∑                          n              =              1                        N                    ⁢                                    w              n                        ⁢                                                            f                  ^                                n                            ⁡                              (                                  x                  ,                  y                                )                                                                        (        13        )            
where N is a number of filters, {circumflex over (f)}n(x,y) is an SSR output image taken through the nth filter (n=1,2, . . . ,N), and wn is a weighting factor. The brightness range of the MSR output image {circumflex over (f)}MSR(x,y) is adjusted to that of the output device by applying the gain/offset correction. FIG. 6 illustrates an MSR output image obtained from a weighted sum of the SSR output images of the input image of FIG. 2A with the LPFs of which support regions are 15, 150, and 500 pixels, respectively. Here, the constants a and b for the gain/offset correction in equation (12) are set to 192 and −30 as in equation 8. As shown in FIG. 6, since the advantageous factors of the SSR output images taken through the filters having different support regions are partially combined, both the global and local contrasts are well increased. Although the halo effect is attenuated to some extent, the MSR output image has a residual halo effect because the most SSR output images have halo effects.
FIG. 7 is a block diagram illustrating the Kimmel method. The Kimmel method is disclosed in the U.S. Pat. No. 6,941,028 and by R. Kimmel, M. Elad, D. Shaked, R. Keshet, and I. Sobel in “A Variational Framework for Retinex,” Int. J. Comput. Vis., vol. 52, no. 1, pp. 7-23, January 2003. In this method, the illumination is iteratively estimated using a QP 701 under the constraints that the illumination varies slowly in all regions including boundaries and its brightness is greater than or equal to the input image etc. After the illumination is estimated, the reflectance is estimated by passing the estimated illumination {circumflex over (l)}(x,y) estimated according to equation (8) through the inverter 707 and multiplying the output of the inverter 707 with the input image f(x,y) at the multiplier 709. The estimated illumination {circumflex over (l)}(x,y) is gamma-corrected at the corrector 703 for decreasing the dynamic range. Finally, the output image is obtained by multiplying the gamma-corrected illumination with the estimated reflectance at the multiplier 705 as following equation 14:{circumflex over (f)}(x,y)={circumflex over (l)}(x,y)α·{circumflex over (r)}(x,y)  (14)
where {circumflex over (r)}(x,y) is an estimated reflectance. This method is advantageous in estimating the reflectance satisfying the constraints. However, the iterative estimation of illumination causes massive computation and no enhancement operation is applied to the {circumflex over (r)}(x,y).
In the color image enhancement methods based on the conventional image formation model, the image enhancement in the RGB color space can be expressed by equation (15) by applying the image of equation (1) to the respective component images fi(x,y),iε{R,G,B}.fi(x,y)=li(x,y)·ri(x,y),iε{R,G,B}  (15)
On the basis of equation (15), the image enhancement can be applied to each component image of the RGB color image. In this case, the red (R), green (G), and blue (B) components vary independently such that the color variations may be caused by variations of the ratio of the R, G, and B components. The color image enhancement methods are proposed in SSR (D. J. Jobson, Z. Rahman, and G. A. Woodell, “Properties and Performance of a Center/surround Retinex,” IEEE Trans. Image Process., vol. 6, no. 3, pp. 451-462, March 1997), MSR, and MSR with Color Restoration (MSRCR) (D. J. Jobson, Z. Rahman, and G. A. Woodell, “A Multi-scale Retinex for Bridging the Gap Between Color Images and the Human Observation of Scenes,” IEEE Trans. Image Process., vol. 6, no. 7, pp. 965-976, July 1997). These methods enable obtaining output color images {circumflex over (f)}i(x,y),iε{R,G,B} only with the enhanced reflectance of equation (16).{circumflex over (f)}i(x,y)=({circumflex over (l)}i(x,y))0·log({circumflex over (r)}i(x,y)),iε{R,G,B}  (16)
In equation (16), it is noted that the output image {circumflex over (f)}i(x,y),iε{R,G,B} is obtained regardless of the estimated illumination of the input image.
In a case that the estimated reflectance of the R, G, and B component images have similar values, the R, G, B component ratio becomes 1:1:1 such that a gray-world violation occurs, in which the color of the output image shifts to gray.
In order to overcome this problem, MSRCR is proposed. In the MSRCR, a color restoration process is added as following equation (17):{circumflex over (f)}MSRCRi(x,y)=Ci(x,y)·{circumflex over (f)}MSRi(x,y),iε{R,G,B}  (17)
where {circumflex over (f)}MSRCRi(x,y) and {circumflex over (f)}MSRi(x,y) are MSRCR output image and an MSR output image of the component image {circumflex over (f)}i(x,y),iε{R,G,B}, and Ci(x,y) is a weighting function expressed as equation 18 for color restoration:
                                                        C              i                        ⁡                          (                              x                ,                y                            )                                =                      c            ·                          log              [                              d                ·                                                                            f                      i                                        ⁡                                          (                                              x                        ,                        y                                            )                                                                                                  ∑                                              j                        =                        1                                            S                                        ⁢                                                                  f                        j                                            ⁡                                              (                                                  x                          ,                          y                                                )                                                                                                        ]                                      ,                  i          ∈                      {                          R              ,              G              ,              B                        }                                              (        18        )            
where c and d are constants for color restoration, S is the number of the color components of the image. In equations 17 and 18, the RGB component images of the MSRCR output color image are weighted according to the RGB component ratio of the input color image. Accordingly, the output image is presented with a color similar to that of the input color image. Finally, the gain/offset correction of equation 12 is applied to the MSRCR output image {circumflex over (f)}i(x,y),iε{R,G,B} to fit the brightness value range of the output image to that of the output device. FIG. 8A illustrates an input color image, and FIGS. 8B and 8C illustrate MSR and MSRCR output color images for the input color image of FIG. 8A. The constants a and b for the gain/offset correction of equation 12 are set to 192 and −30, respectively, and the constants c and d for the color restoration of equation 18 are set to 46 and 125, respectively, for the MSRCR as proposed by D. J. Jobson, Z. Rahman, and G. A. Woodell, in “A Multi-scale Retinex for Bridging the Gap Between Color Images and the Human Observation of Scenes,” IEEE Trans. Image Process., vol. 6, no. 7, pp. 965-976, July 1997. In the MSR output color image of FIG. 8B, it is shown that the contrast of the dark area in the input color image is well enhanced, but the gray-world violation appears. Although the gray-world violation is attenuated to an extent in the MSRCR output color image of FIG. 8C, residual gray-world violation appears. That is, the MSRCR cannot be a solution for avoiding the gray-world violation effect occurred in the MSR.